
*********************
* BOUNDING ANALYSIS *
*********************
* Labussiere 2023

/* Produce the following outputs:
	In main text:
		Table 4
	
	In supplementary materials:
		Table S12
*/


use full_sample.dta


* Rename variables for clarity 
gen outcome = cito_score_st_year

gen age_nat_dummies = age_nat_13_censored
label values age_nat_dummies agenat_15

gen age_nat_cat = age_nat_13_cat7
label values age_nat_cat agecat_7


capture log close
log using psacalc_1000.smcl, replace 

* sibling fixed-effects (-> to get rho)
xtreg outcome i.nat_13_censored i.first_born_TOT i.gender ///
ib9.date_birth_cat i.homeowner_alo_cito i.q3_st_disp_income_cito ///
i.secm_ma_agg2_cito i.secm_pa_agg2_cito b(1).type_household_cito ///
nbr_children_cito ///
if ind_in_SFE, fe

scalar rho = e(rho)
	* 0.5128
	
* full OLS model (-> to get R^2)
regress outcome i.nat_13_censored i.gender i.first_born_TOT i.date_birth_cat ///
	i.parent_ed_level_miss i.homeowner_alo_cito i.q3_st_disp_income_cito ///
	i.secm_ma_agg2_cito i.secm_pa_agg2_cito i.country_of_birth_ma_agg ///
	ib(1).type_household_cito nbr_children_cito, /*coeflegend*/
estimates store OLS
	
scalar O_R2max = 1.3*e(r2)
	* 0.1274

* Macro for mcontrol() option
global def_mcontrol mcontrol(1.gender 1.first_born_TOT ///
						2.date_birth_cat 3.date_birth_cat 4.date_birth_cat ///
						5.date_birth_cat 6.date_birth_cat 7.date_birth_cat ///
						8.date_birth_cat 9.date_birth_cat)
	

******************************************************************
* TABLE 4
* Bounding analysis: calculating Beta with Rmax = 1.3 R2 and  
* varying values for Delta
******************************************************************

capture program drop varying_delta_O_R2max
program define varying_delta_O_R2max
bs r(beta), rep(1000): psacalc beta 1.nat_13_censored, delta(`1') rmax(`=O_R2max') ///
 model(regress outcome i.nat_13_censored i.gender i.first_born_TOT i.date_birth_cat ///
	i.parent_ed_level_miss i.homeowner_alo_cito i.q3_st_disp_income_cito ///
	i.secm_ma_agg2_cito i.secm_pa_agg2_cito i.country_of_birth_ma_agg ///
	ib(1).type_household_cito nbr_children_cito) `2'
end


varying_delta_O_R2max 0.5 "${def_mcontrol}" 
varying_delta_O_R2max 1 "${def_mcontrol}" 
varying_delta_O_R2max 1.5 "${def_mcontrol}" 
varying_delta_O_R2max 2 "${def_mcontrol}" 
varying_delta_O_R2max 3 "${def_mcontrol}" 



******************************************************************
* TABLE S12
* Bounding analysis with Delta = 1 and varing values for Rmax 
******************************************************************

capture program drop varying_Rmax_1
program define varying_Rmax_1
bs r(beta), rep(1000): psacalc beta 1.nat_13_censored, rmax(`1') ///
 model(regress outcome i.nat_13_censored i.gender i.first_born_TOT i.date_birth_cat ///
	i.parent_ed_level_miss i.homeowner_alo_cito i.q3_st_disp_income_cito ///
	i.secm_ma_agg2_cito i.secm_pa_agg2_cito i.country_of_birth_ma_agg ///
	ib(1).type_household_cito nbr_children_cito) `2'

end

varying_Rmax_1 0.13 "${def_mcontrol}"
varying_Rmax_1 0.15 "${def_mcontrol}"
varying_Rmax_1 0.17 "${def_mcontrol}"
varying_Rmax_1 0.19 "${def_mcontrol}"
varying_Rmax_1 0.2 "${def_mcontrol}"


log close



	
